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A class of strategies has been proposed to reduce contrail formation in the United 
States airspace. A 3D grid based on weather data is built and the cruising altitude level 
of aircraft is adjusted to avoid the persistent contrail potential area with the consideration 
to fuel-efficiency. In this paper, the authors introduce a contrail avoidance strategy on 
3D grid by considering additional operationally feasible constraints from an air traffic 
controliers aspect. First, shifting too many aircraft to the same cruising level will make 
the miles-in-trail at this level smaller than the safety separation threshold. Furthermore, 
the high density of aircraft at one cruising level may exceed the manageable workload 
for the traffic controiler. Therefore, in our new modei we restrict the number of total 
aircraft at each level. Second, the aircraft count variation for successive intervals can 
not be too drastic since the workload to manage climbing/descending aircraft is much 
larger than managing cruising aircraft. The contrail reduction problem is formulated as 
integer programming and the problem is shown to have the property of total unimodularity. 
Solving the corresponding relaxed linear programming with the simplex method provides an 
optimal and integral solution to the problem. Simulation results are provided to illustrate 
the methodology. 


Nomenclature 


T temperature in Celsius 

RHw relative humidity with respect to water 
RHi relative humidity with respect to ice 
l cruising altitude level index 

R\ potential persistent contrail formation matrix (contrail matrix) at time t at level l 
A\ aircraft position matrix at time t at level l 
CFI center contrail frequency index 
A^ f e level changing constraint 

bi(l) the number of aircraft which cruise at level l in the original flight plans 
62(0 center level capacity at level l 

A q given threshold of climbing/descending aircraft count 


I. Introduction 

Contrails are artificial clouds that are the visible trails of condensed water vapor made by the exhaust 
of aircraft engines. Depending on atmospheric conditions, contrails may be visible for only a few seconds or 
minutes, or may persist for many hours which may affect climate. Persistent contrails reduce incoming solar 
radiation and outgoing thermal radiation in a way that accumulates heat. The global mean contrail cover 
in 1992 was estimated to double by 2015, and quadruple by 2050 due to air traffic increase. Studies suggest 
that the environmental impact from persistent contrail is estimated to be three or four times, or even ten 

* School of Aeronautics and Astronautics, Purdue University 
^NASA Ames Research Center, Moffett Field, CA 


1 of 13 


Copyright © 2012 by Peng Wei. Published by the American Institute of Aeronautics and Astronautics, Inc., with permission. 


times larger than aviation-induced emissions. Therefore, strategies and policies to reduce aircraft induced 
persistent contrail need to be studied to minimize the impact on the global environment. 

There have been some research works accomplished to identify and reduce persistent contrail forma- 
tion. Gierens and Noppel 1 reviewed various strategies for contrail avoidance including changing engine 
architecture, enhancing airframe and engine integration, using alternate fuels, and modifying traffic flow 
management procedures. Among the traffic flow management solutions, Mannstein presented a strategy to 
reduce the environment impact of contrails significantly by small changes to each aircraft’s flight altitude. 
Campbell ' proposed a mixed integer programming method to optimally reroute aircraft to avoid the forma- 
tion of persistent contrails. Both methods require the onboard contrail detection system. Fichter showed 
that the global annual mean contrail coverage could be reduced by decreasing the aircraft cruise altitude. 
Williams ’ proposed strategies for contrail reduction by identifying fixed and varying maximum altitude 
restrictions. These restrictions generally require more fuel burn and add congestion to the already crowded 
airspace at lower altitudes. Sridhar et al. provided a set of strategies to reduce contrail formation. A 3D 
grid model was constructed based on weather data to describe the contrail potential areas and the cruising 
level of aircraft is adjusted to avoid these areas with the consideration of fuel-efficiency. 

The goal of this paper is to refine the research work in by introducing two operationally feasible 
constraints from air traffic control perspective with the consultation from experienced and expert controllers. 
Moreover, we formulate the problem as integer programming. Because of the total unimodularity of its 
constraint matrix, solving the relaxed linear programming by simplex method will yield the optimal and 
integral solution with much shorter computation time. 

The rest of the paper is organized as follows. Section II introduces the atmospheric data, contrail 
formation model, aircraft data and contrail frequency index used in this paper. Section III presents the 
problem formulation and our contrail reduction strategy. The numerical results are shown in Section IV. 
Section V summarizes and concludes the paper. 


II. Data and Model 

A. Atmospheric Data 
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Figure 1. Contours of temperature and RHw at 34,057 feet at 8AM EDT on August 1, 2007. 

Contrails can be observed from surface data and detected by satellite data. Duda 1 ' has related 
the observations to numerical weather analysis output and showed that persistent contrail formation can be 
computed using atmospheric temperature and humidity data retrieved from the Rapid Updated Cycle (RUC) 
data provided by the National Oceanic and Atmospheric Administration (NOAA). Contrails can persist when 
ambient air is supersaturated with respect to ice, which means that relative humidity with respect to ice 
(RHi) is greater than one hundred percent. The RHi can be computed from relative humidity with respect 
to water (RHw) and temperature, which are available in the RUC data. The one-hour forecast and the 
40-km RUC data are used in this paper. The data have a temporal resolution of one hour, a horizontal 
resolution of 40 kilometers, and isobaric pressure levels from 100 to 1000 hectopascals (hPa) with 25 hPa 
increments. The vertical range of interest in this study is from 150 hPa to 400 hPa, which is equivalent to 
pressure altitude of about 23,600 feet to 44,400 feet. The temperature and RHw contours at 8AM eastern 
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daylight time (EDT) on August 1, 2007 at pressure altitude 250 hPa, or 34,057 feet, are shown as the left 
and right subfigures in Fig. 1. 


B. Contrail Formation Model 

The potential persistent contrail formation areas (contrail areas) are defined as areas with RHi greater than 
or equal to 100%. RHi can be computed from RHw and temperature using the saturation vapor pressure 
coefficients of Alduchov, formulated as 


RHi = RHw x 


6.0612e 181O2T/(249 - 52+T) 
6.1162e 22 - 577T /( 237 - 78+T ) ’ 


( 1 ) 


where T is the temperature in Celsius. The atmospheric profile shown in the left and right subfigures in 
Fig. 1 can be translated to a contour of RHi, as shown in Fig. 2. 
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Figure 2. Contour of RHi at 34,057 feet at 8AM EDT on August 1, 2007. 


The 40-km RUC data have (113 x 151) data points. The altitude level index l is defined as l = 1, 2, ..., 11 
corresponding to isobaric pressure level at 400,375, ..., 150 hPa. The level index, isobaric pressure level, and 
approximate aircraft cruising altitude are listed in Table 1. 


Table 1. Level index, isobaric pressure level and approximate aircraft cruising altitude. 


Level index 

1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

11 

Pressure level (hPa) 

400 

375 

350 

325 

300 

275 

250 

225 

200 

175 

150 

Cruising altitude (100 feet) 

236 

251 

267 

283 

301 

320 

341 

363 

387 

414 

444 


The potential persistent contrail formation matrix (contrail matrix) at time t at level l is defined as 


( r i,i 


R't = 


r 1,2 


?T,151 


\ 


( 2 ) 


\ 113,1 1* 113,2 r 113,151 / 


where r lt j is 1 if RHi > 100% at grid (i,j) and 0 if RHi < 100%. 


C. Aircraft Data 

Contrails form when aircraft fly through a potential contrail formation area. Thus aircraft locations are 
needed to determine the contrail formation frequency. The aircraft data used in this paper are obtained 
from the aircraft locations provided by the Federal Aviation Administration’s (FAA’s) Aircraft Situation 
Display to Industry (ASDI) data. The ASDI has a sampling rate of one minute. The same 3D grid used in 
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the RUC data is used to generate the aircraft position matrix. The aircraft position matrix is defined as 


K = 


Ol,l 


ai,2 


Ol,151 


\ 0113, 1 0113,2 


( 3 ) 


0113,151 J 


where dij is the number of aircraft within grid (i,j) flying closest to level l at time t. The aircraft position 
matrix indicates the air traffic density in the grid scale at different altitudes. 


D. Contrail Frequency Index 

Contrail frequency index is defined as the number of aircraft that would fly through potential contrail 
formation regions during a period of time. Center contrail frequency index is used to indicate the contrail 
severity in a given center. s To specify the location of the twenty U.S. centers in the grid scale, the center 
grid matrix is defined as 


( 


C a 


Cl,l C i,2 


Cl, 151 


\ 


( 4 ) 


\ Cn3,i Cn3,2 • ■ ■ Cn3,i5i J 


where Cj,j is one if the grid point is within the center and zero if not. 

The center contrail frequency index is defined as the number of aircraft flying through contrail area at 
time t at level l, formulated as 

113 151 

C F I center, l,t — (S) 

i = 1 1 = 1 

where Tij, dij, Cij are defined in Eqn. (2), (3), (4). 

The contrail frequency index derived in the previous section indicates the actual contrail activities. For 
strategic planning, prediction of the contrail frequency at time t at level l of a certain center is calculated 

by 

113 151 

C F ' I center, l, t = («> 

*= i i=i 

where r( 7 is from RUC forecast data and a[ j is the predicted aircraft locations from historical air traffic 
data. 


III. Problem Formulation and Contrail Reduction Strategy 

The center contrail frequency index can be used to identify the flight level that would have formed the 
most contrails and find an alternate cruising altitude with less contrail formations. The contrail frequency 
index after the contrail reduction strategy has been applied is formulated as 

113 151 

CF I center, l,t = EEw* in 

*= i i=i 

where r» j and c; j are defined in Eqn. (2) and (4) and di j is defined in Eqn. (3) with the aircraft location after 
applying the contrail reduction strategy, which can be considered as variables. Thus the contrail reduction 
strategy is to solve the optimization problem, whose objective is to minimize CFI under several constraints. 
In this paper, the authors consider the fuel efficiency constraint and the operationally feasible constraints. 

A. Fuel Efficiency Constraint 

As discussed in, 12 the original cruising altitude filed in flight plan is usually the optimal cruising altitude in 
term of fuel efficienccy and aircraft performance. Therefore shifting aircraft up and down to other further 
levels is more fuel consuming than flying them at original altitude. In this paper, we allow aircraft to 
alternate the cruise level within the range (block altitude clearance 1 ' 1 ) [l — A lf e ,l + A lf e \, where l is the 
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original cruising level and A l f e is the level changing constraint based on fuel efficiency. A l f e is set to 1 and 
2 respectively in our simulation. 

As an example, the contrail freqency index matrix at Atlanta Center at 8AM EDT on August 1, 2007 is 


CFL 


center, t — 


0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

1 

0 

0 

0 

0 

2 

0 

0 

1 

0 

28 

26 

31 

33 

10 

52 

68 

105 

51 

6 

7 

57 

74 

43 

25 

13 

76 

104 

148 

71 

5 

10 

69 

29 

44 

33 

28 

105 

128 

209 

132 

9 

14 

16 

12 

5 

8 

5 

62 

47 

36 

22 

6 

0 

0 

2 

0 

0 

5 

65 

45 

35 

19 

6 

0 

0 

0 

0 

0 

4 

71 

108 

36 

19 

6 

0 

0 

0 

0 

0 

10 

74 

122 

33 

19 

0 

0 


( 8 ) 


where the bold diagonal item CFI l centert ; is the contrail freqency index at level l before applying contrail 
reduction and the off diagonal item CFI l center t ; at (/', l ) is the contrail freqency index when guiding aircraft 
from level l to level l'. There are 11 levels in total as illustrated in Table 1, therefore the matrix (8) is 11 x 11. 


only allowed to move 

two levels 

up or 

down (A l f e 

= 2) tl 


/ 0 

0 

0 

X 

X 

X 

X 

X 

X 

X 

X 


0 

0 

0 

0 

X 

X 

X 

X 

X 

X 

X 


0 

0 

0 

0 

0 

X 

X 

X 

X 

X 

X 


X 

1 

0 

0 

0 

0 

X 

X 

X 

X 

X 


X 

X 

31 

33 

10 

52 

68 

X 

X 

X 

X 


X 

X 

X 

25 

13 

76 

104 

148 

X 

X 

X 


X 

X 

X 

X 

28 

105 

128 

209 

132 

X 

X 


X 

X 

X 

X 

X 

62 

47 

36 

22 

6 

X 


X 

X 

X 

X 

X 

X 

45 

35 

19 

6 

0 


X 

X 

X 

X 

X 

X 

X 

36 

19 

6 

0 


V X 

X 

X 

X 

X 

X 

X 

X 

19 

0 

0 

invalid 

move 

More strictly, if the aircraft 

can only be 


(9) 


(A lf e = 1) the matrix is 


0 

0 

X 

X 

X 

X 

X 

X 

X 

X 

X 

0 

0 

0 

X 

X 

X 

X 

X 

X 

X 

X 

X 

0 

0 

0 

X 

X 

X 

X 

X 

X 

X 

X 

X 

0 

0 

0 

X 

X 

X 

X 

X 

X 

X 

X 

X 

33 

10 

52 

X 

X 

X 

X 

X 

X 

X 

X 

X 

13 

76 

104 

X 

X 

X 

X 

X 

X 

X 

X 

X 

105 

128 

209 

X 

X 

X 

X 

X 

X 

X 

X 

X 

47 

36 

22 

X 

X 

X 

X 

X 

X 

X 

X 

X 

35 

19 

6 

X 

X 

X 

X 

X 

X 

X 

X 

X 

19 

6 

0 

X 

X 

X 

X 

X 

X 

X 

X 

X 

0 

0 


x \ 


(10) 


The following two constraints are called operationally feasible constraints for contrail reduction. 
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B. Center Level Capacity Constraint 

In, 1 the authors shifted all the aircraft at level l to the level with minimum CFI l centert v However, if there 
is one level always with the minimum CFI l center f ; , the previous approach may potentially form a congested 
center level by moving a large number of aircraft to the same level. In practice, the wake turbulence or wake 
vortex phenomenon constrains the air traffic controller to maintain a minimum miles-in-trail. The wake- 
turbulence separation criterion is currently a limiting factor in airspace capacities. The FAA is working with 
NASA to develop and demonstrate integrated systems technology for addressing separation criteria. Some 
work has been done on miles-in-trail separation and space metering. ' ' In summary, it is potentially 

harmful to move too many aircraft to the same level because of smaller separations and associated higher 
risk. In Fig. 3, the miles-in-trail (MIT) changes from MIT1 to a smaller value MIT3 after the metering point 
because there are some aircraft merging from the lower level. 

>K > • >)- >)- 

Metering Point 

X 

X MIT2 

Figure 3. Miles-in-trail is reduced after more aircraft merge from the lower level. 


Another reason to introduce this constraint is to keep the air traffic controller’s workload below a reason- 
able threshold. In order to measure dynamic density and evaluate the controller’s workload, several metrics 
were proposed. J ’ Among all these dynamic density metrics, the aircraft count or the traffic density in one 
area is the first factor noticed by every metric. Therefore maintaining the aircraft count at each level below 
a threshold value will guarantee the workload of air traffic controller manageable. 

In this paper, we propose the center level capacity constraint to keep the aircraft amount at each level 
lower than a critical value. This idea is similar to the Monitor Alert Parameter of a sector, which is set 
to reflect controller’s acceptable workload. We use our level capacity constraint as a straightforward and 
operationally intuitive method to keep the controller workload at a manageable level. 

C. Climbing/Descending Aircraft Count Constraint 



Figure 4. Aircraft are climbing from cruising level l to level l + 1. 


Normally it is difficult for the controllers to visualize potential vertical conflicts 


with a 2-D radar scope. 
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As shown in Fig. 4, controllers can observe and manage the horizontal aircraft (in black) at level l and level 
l + 1 respectively. However, if there are too many climbing aircraft (in blue) moving from level l to level 
l + 1 , it will introduce extra workload to air traffic controllers. 

Almost all the metrics in • mentioned the number of climbing/descending aircraft or the number of 
aircraft with large altitude change as the dynamic density impact factor, which means a huge number of 
climbing/descending will increase the controllers’ workload. 

In this study, suppose the solution of aircraft count at level l at time t is q* enter fi , we restrict the 
differences of current solution with the aircraft counts at time t — 1 and t + 1 within a given threshold, i.e. , 
k*center,t,i ~ q C enter,t-i,i\ < A q and | q* cen ter,t,i “ <?center,t+i,i | < Ag. By doing this, our contrail reduction 
strategy should not shift too many aircraft up and down compared to the original flight plans. Therefore 
less extra workload will be added to the controllers. 


D. Problem Formulation 

Instead of choosing the aircraft count in each grid as the variables, we use the aggregate aircraft count at 
each level to take advantage of the special structure of contrail frequency index matrix in Eqn. ( 8 ). 

For ease of illustration, A lf e is set to 1 at this time. Thus our variables are chosen based on the +/ — 1 
level changing scheme in Eqn. (10), and the unknown variables are arranged as a vector of length 31 

[xi,l,X2,l,Xl,2, X2,2,X3,2, •••! ^9,10; ^10,10 j ^11,10) ^10,11) 2ftl, Ilf; 


where Xij is the number of aircraft which used to cruise at level j and now cruise at level i. 

Our goal is to guide aircraft to the level with minimum CFI value while satisfying all the constraints. 
According to Eqn. (10), the c vector in the objective function is: 


c = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 33, 0, 10, 13, 52, 76, 105, 104, 128, 47, 209, 36, 35, 22, 19, 19, 6, 6, 0, 0, 0]'. 


Fuel Efficiency Constraint tells us the aircraft used to cruise at level j can only be adjusted to level 
j — 1, j and j + 1 (except the aircraft used to cruise at the lowest level 1 and the highest level 11). By flow 
conservation, we have matrix A\ of 11 x 31 and vector bi of length 11 : 


A i 


/ 1 1 0 0 0 
0 0 111 

0 0 0 0 0 

V 0 0 0 0 0 


0 0 0 0 0 \ 

0 0 0 0 0 

1110 0 

0 0 0 1 1 / 


bi 


[10, 10, 10, 15, 40, 200, 200, 100, 50, 20, 10]', 


( 11 ) 


where bi(l) is the number of aircraft which used to cruise at level l in the original flight plans. 

Center Level Capacity Constraint requires the total aircraft number at each level after contrail reduction 
is fewer than the capacity. Since center level capacity is not an implemented setting in practice, we estimated 
the values of the capacities for each level. For example, there can be three source streams of aircraft forming 
the aircraft cruising at level l , which are x^i - 1 aircraft from level l — 1 , Xi t i aircraft staying at level l and 
Xi : i+i aircraft shifted from level l + 1. In summary, we obtain constraint matrix A 2 of 11 x 31 and vector 62 
of length 11 : 

/ 1 0 1 0 0 0 ... 0 0 0 0 0 0 \ 

010101 ... 000000 


A 2 = 


( 12 ) 


0 0 0 0 0 0 

V 0 0 0 0 0 0 


10 10 10 
0 0 0 1 0 1 / 


b 2 = [ 15 , 15, 20, 25, 50, 250, 250, 150, 


70,30,15]', 


where b 2 (l) is the center level capacity at level l. 
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Climbing/Descending Aircraft Count Constraint requires the absolute difference between the adjusted 
aircraft at level l and the aircraft used to cruise at level l at time t — 1 smaller than a given threshold, i.e., 
\q* enter, t,i ~ Qcenter,t-i,i\ < A q. Similarly, the absolute difference between q* centert l and q ce nter,t+i,i should 
also be less or equal than A q. 

It is easy to calculate the number of aircraft used to cruise at level l at time t—1 and t + 1 from ASDI 
aircraft data. However, A q value is not necessarily the same under different altitudes, center or sector 
boundary shapes and geographic conditions. In the simulation part of this paper we assume all the A q are 
the same. 

Ignoring the center subscript notion, we denote the number of aircraft used to cruise at t — 1 and t + 1 
at level l as qt-i,i and qt+i i, which we can compute from the historical data. At each level, we restrict the 
maximal number of climbing/descending aircraft with A qi. Then we have: 


-Aqi 

< 

*1,1 + *1,2 - qt+ 1,1 

< 

Agi, 

-A q 2 

< 

*2,1 + *2,2 + *2,3 — qt+ 1,2 

< 

Ag 2 , 

-Agn 

< 

*11,10 + *11,11 — <Zt+l,ll 

< 

Agn, 

-Agi 

< 

*1,1 + *1,2 — qt- 1,1 

< 

Agi, 

-A q 2 

< 

*2,1 +*2,2 +*2,3 ~qt- 1,2 

< 

A g 2 , 

-Agn 

< 

*11,10 + *11,11 — <7t-l,ll 

< 

Agn. 


Let’s denote A 2 as A t+ i i+ , (q t+1 + A q) as b t + 1 ,+ , -A 2 as A t+1 _, (~q t + 1 + A q) as b t+1 _, A 2 as A t - i >+ , 
( qt-i + A q) as b t - 1 ,+, — A 2 as A t - 1 ,_, (— q t +i + A q) as b t -i~, the Climbing/Descending Aircraft Count 
Constraint can be written as: 


( A t+ 1,+ \ 


( &*+!,+ \ 

At+i- 

X < 

bt+ 1,_ 

At~ i.+ 


b t - 1,+ 

V ^t-b- } 


l h- 1,- / 


By combining (11) (12) (13) , we have the fule efficiency constraint (one equality) and two operationally 
feasible constraints (five inequalities in total) listed all as inequalities below: 


l A! \ 


f bl 

—Ai 


-bi 

a 2 


b 2 

At+ 1,+ 

X < 

bt+ 1,+ 

At+i,_ 


b t + 1- 

At- 1,+ 


b t - 1,+ 

\ A t - 1,_ / 


V bt- 1,- / 


which is simply denoted as: 


Ax < b. 

Finally, we formulate an integer programming problem: 

min c' x 


s.t. Ax < b, 

x > 0, and x £ L 


(14) 


(15) 
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E. Solution — Refined Contrail Reduction Strategy 

In this work the integer programming (15) is relaxed to a linear programming (16) for computational effi- 
ciency. 


min c'x 

s.t. Ax < b, (16) 

x > 0 

We study the total unimodularity of matrix A in (16) and prove that there exists an optimal and integral 
solution for the linear programming relaxation, which is also the optimal solution for the original integer 
programming (15). It is guaranteed to be optimal and integral when (16) is solved by the simplex method. 

Theorem 1 Matrix A in (16) is total unimodular. 

Proof 1 Here the dimension of matrix A is 66 x 31. Targeting to transform each row of submatrix A 2 (14) 
into the row with a single 1, through a series of elementary column operations, ' we have transformed A as: 


1 

1 

-1 

-1 

0 

-1 

0 

0 ... 

0 

0 1 

0 

0 

1 

1 

1 

0 

-1 

0 ... 

0 

0 

0 

0 

0 

0 

0 

1 

1 

1 ... 

0 

0 

o • • 

0 

0 

0 

0 

0 

0 

0 ... 

-1 

0 

0 

0 

0 

0 

0 

0 

0 

0 ... 

0 

-1 

0 

0 

0 

0 

0 

0 

0 

0 ... 

1 

1 

-1 

-1 

1 

1 

0 

1 

0 

0 ... 

0 

0 

0 

0 

-1 

-1 

-1 

0 

1 

0 ... 

0 

0 

• o 

0 

0 

0 

0 

-1 

-1 

-1 ... 

0 

0 

o • • 

0 

0 

0 

0 

0 

0 

0 ... 

1 

0 

0 

0 

0 

0 

0 

0 

0 

0 ... 

0 

1 

0 

0 

0 

0 

0 

0 

0 

0 ... 

-1 

-1 

1 

0 

0 

0 

0 

0 

0 

0 ... 

0 

0 

0 

1 

0 

0 

0 

0 

0 

0 ... 

0 

0 

• o 

0 

0 

0 

1 

0 

0 

0 ... 

0 

0 

o • • 

0 

0 

0 

0 

0 

0 

0 ... 

1 

0 












\ ■■■■■■■■■■■ 


where the part between the second and the third horizontal lines is the transformed submatrix A 2 . Since 
matrices A t +i t +, A t+ 1 : _ ; A t - i.+, A t _ 1 ,_ in (14) are either A 2 or —A 2 , every line in the omitted part below 
the third horizontal line also only contains one single 1 or one single —1. 

Based on (f3)(v) of, if the transformed submatrices Ai and —A\ above the second horizontal line in 
(17) is totally unimodular, then the original matrix A is also totally unimodidar. Moreover, if the transformed 
Ai is total unimodular, then the two parts of the transformed A 1 and —A\ above the second horizontal line 
are total unimodular. 
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Now let's prove the following matrix is total unimodular: 


0 ... 0 

0 ... 0 

1 ... 0 


/ 1 1 —1 —1 0 — 1 0 

001 110-1 
0 0 0 0 0 1 1 


0 ^ 

0 

0 


0 0 
0 0 
\0 0 


0 0 
0 0 
0 0 


0 0 
0 0 
0 0 


0 0 
0 0 
0 0 


-1 

0 

1 


0 

-1 

1 


(18) 


where each column either contains only one number 1 or a 1/ — 1 pair. By (43) (v) of, the first two columns 
with one non-zero item 1 can be removed while preserving the total unimodularity of a matrix. Then according 
to (18) of, the matrix whose each column contains exactly one 1 and exactly one 1 is totally unimodular. 


Since matrix A is totally unimodular and vector b is integral, there must exist an optimal and integral 
solution for both the linear programming relaxation (16) and for the integer programming in (15), which is 
guaranteed by any simplex method. 

In this paper we are dealing with a medium size linear programming problem, so the simplex method 
in Matlab is capable to provide the efficient solution. The detailed simulation results are presented in next 
section. 


IV. Simulation Results 

In this section, we use the same RUC weather data and ASDI aircraft data as in (8), (9), (10). We first 
evaluate the performance of our refined strategy after introducing addtional operationally feasible constraints. 
Secondly, the aircraft count at each cruising level is examined to see if the result of our refined strategy can 
be governed by the estimated level capacity. Thirdly the numbers of aircraft shifted to each level are plotted 
which can also be understood as climbing/descending aircraft counts. 

A. Contrail Reduction Performance 

In Fig. 5, the indice 1, 2, 3 on x-axis represent the historical data (no contrail reduction), the original contrail 
reduction strategy in and our refined strategy. The red curve is the result of contrail reduction by allowing 
+/ — 1 level aircraft shifting. The blue curve is the contrail reduction by allowing +/ — 2 levels aircraft 
shifting. 
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Figure 5. Contrail reduction performance comparison. 
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The original strategy has better contrail reduction performance whenever it is the +/ — 1 level shifting 
or the +/ — 2 levels shifting. The +/ — 2 levels shifting of the original strategy gives out almost 60% contrail 
reduction while the +/ — 2 levels shifting of our strategy can provide a 20% contrail reduction. 

B. Level Capacity Constraint 

In Fig. 6, the numbers 1, 2, ..., 11 on x-axis are the 11 cruising levels and the j/-axis indicates the number 
of aircraft cruising at each level. The red curve describes the level capacity which helps the controllers keep 
their workload manageable. The blue curve is the adjusted aircraft count at each level provided by the 
original strategy. The purple curve is the contrail reduction result provided by our refined strategy. The left 
subfigure is the result of the +/ — 1 level shifting while the right subfigure is the result of the +/ — 2 levels 
shifting. 




Figure 6. Our refined strategy is strictly governed by the level capacities. 


The contrail reduction result from the original strategy will potentially guide the aircraft to the level 
with minimal contrail potential grids, which may form a very busy cruising level, i.e., level 9 in the right 
subfigure of Fig. 6. However, our strategy successfully adjust the aircraft to reduce the contrail frequency 
index and restrict the aircraft count below the level capacity at the same time. 

C. Shifted Aircraft Number 
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Figure 7. The climbing/descending aircraft count at each level. 


Fig. 7 shows us the shifted aircraft number at each level, which can also be understood as climb- 
ing/descending aircraft counts. The blue curve represent the result of the original strategy and the purple 
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curve is provided by our refined strategy. The left subfigure is the result of the +/ — 1 level shifting while 
the right subfigure is the result of the +/ — 2 levels shifting. 

No matter it is the +/ — 1 level shifting or the +/ — 2 levels shifting, our strategy does great job on 
managing the aircraft count of each level at a relative low value. This feature can help the controllers process 
fewer climbing/descending aircraft at each level. 

V. Conclusion 

This study brings two additional operationally feasible constraints to the original contrail reduction 
problem. By introducing the center level capacity constraint and the climbing/descending aircraft count 
constraint, the paper formulates the 3D grid contrail reduction problem as an integer programming. With 
the total unimodular property, the relaxed linear programming provides the optimal and integral solution in 
a computational efficient manner, which is the most important advantage of our approach. Another benefit 
of the method is that instead of shifting all the aircraft from one level to another, our formulation supports 
the feature of splitting the total number of aircraft used to cruise at level l to multiple levels, i.e. , part 
of the aircraft used to cruise at level l can stay at this level, part of them may be shifted to upper levels, 
and the remaining part will be guided to lower allowed cruising levels. Simulation results show that the 
refined contrail reduction strategy provides less contrail reduction than the original strategy presented in. 
However, it strictly obeys the level capacity and the numbers of shifted aircraft at each level are restricted 
to maintain air traffic controllers workload. Future work will integrate operationally acceptable rerouting 
technique ' ’ into contrail avoidance trajectory optimization. Instead of the 3D grid aircraft shifting, the 
new avoidance trajectories and flight plans may provide more contrail reductions. 
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